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Abstract. The phonon density of states (DOS) of graphene with different types of 
point defects (carbon isotopes, substitution atoms, vacancies) is considered. Using a 
solvable model which is based on the harmonic approximation and the assumption that 
the elastic forces act only between nearest neighboring ions we calculate corrections 
to graphene DOS dependent on type and concentration of defects. In particular 
the correction due to isotopic dimers is determined. It is shown that a relatively 
small concentration of defects may lead to significant and specific changes in the 
DOS, especially at low frequencies, near the Van Hove points and in the vicinity of 
the K-points of the Brillouin zone. In some cases defects generate one or several 
narrow gaps near the critical points of the phonon DOS as well as resonance states 
in the Brillouin zone regular points. All types of defects are characterized by the 
appearance of one or more additional Van Hove peaks near the (Dirac) K points 
and their singular contribution may be comparable with the effect of electron-phonon 
interaction. Besides, for low frequencies and near the critical points the relative change 
in density of states may be many times higher than the concentration of defects. 
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1. Introduction 

Due to the low atomic mass of carbon and high binding energy of the valent sp2 bonds 
PQ, graphene-based nanostructures (single- or several-layered graphene and nanotubes 
of different kinds) have high rigidity in one direction combined with excellent flexibility 
in the others [21 [3] leading to the extraordinary sound velocity (about 20000 km/sec) 
[UEIE] and thermal conductivity [H [8j [9j [10] . Due to these properties carbon compounds 
can be used not only for plastic and hydrocarbon resin reinforcement, but also as one 
of the basis materials for nano-mechanics and nano-electronics [111 E2] ■ 

It is natural to expect that as well as in bulk crystals even small concentrations of 
point defects in graphene-based ID and 2D nanostructures (directly observed in [131 E]) 
may lead to specific shifts, broadenings and additional characteristic singularities in 
the electron and phonon densities of states and thus change their optical absorption, 
low-temperature specific heat and transport properties. Such effects in graphene and 
carbon nanotubes are significant because of the occurrence of three isotopes in natural 
carbon ( 12 C, 13 C, 14 C) with the part of the 13 C isotope in the chemically pure carbon 
exceeding one per cent and also because of the high solubility of substitution defects 
of trivalent atoms (such as aluminium, boron and nitrogen) and monovalent atom 
adsorption susceptibility According to [15] , the defect density in a graphene monolayer 
stabilized on a substrate can reach several per cent and for some applications it can be 
additionally doped to raise the electrical conductivity. It is shown that unintentional 
doping of pristine unprocessed graphene under ambient conditions can reach as high 
as 1% [16], while the highest achieved doping level of N is about 5% [Hj. It was also 
computationally established that the graphene planar structure is kept even for 12% Al 
and 20% N concentrations [T8l IT9] . In addition, most of the chemically adsorbed atoms 
(especially monovalent ones) can be treated as isotopic defects because they are bound 
to carbon atoms by the 7r-electron bonds which are not involved in lattice formation 
(this may need slight correction of a electrons' binding energies and angles, but as the 
first approximation they can be taken as in the ideal graphene). It is worth mentioning 
that the impact of defects on the electronic properties of graphite, graphene and carbon 
nanotubes along with a detailed investigation of how electron-phonon interaction affects 
the phonon dispersion curves especially near the K points of the Brillouin zone were 
thoroughly studied in a large number of works [HJ El UOl [20l EH [221 E31 [2^] , while too 
little attention has been paid so far to description of the direct influence of defects on 
the phonon spectra. However the anomalies of the phonon spectra due to the electron- 
phonon interaction may be visibly distorted by the Van Hove spikes induced by defects. 

This paper is devoted to the description of the effect of some point defects (isotopic 
defects, substitutional atoms, vacant lattice sites) on the phonon spectrum of graphene. 
As a starting point we consider in Section 2 the ideal graphene phonon spectrum in the 
simple harmonic approximation assuming that elastic forces act only between nearest 
neighboring hard ions and are described by three harmonic force constants J\, J2, J3 
corresponding to three different parts of the interatomic interaction: the central(l) and 
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non-central(2) in-plane forces and the empirical to-the-plane backmoving non-central 
force (3). The values of these constants are chosen to get the least discrepancy of 
the eigenfrequencies calculated in the framework of a simple model for the T, K and 
M-points of the Brillouin zone with those obtained in well-performed inelastic x-ray 
scattering experiments jl]. It is important to underscore that a convergence within the 
limits of experimental error of theoretical and experimental phonon dispersion curves for 
graphene is unlikely to be attainable on the base of the simple three-parameter nearest- 
neighbors model. For example, by using of such a model it is impossible to explain the 
phenomenon of " overbending" , which was observed on the dispersion curves for graphite 
and graphene jH [221 US\ • However, being operationally rather simple, it gives for 
appropriate values of the constants J\, J 2 , J3 the dispersion curves and phonon density 
of states (DOS) for graphene, which are qualitatively just similar and quantitatively 
are rather close to those obtained using theoretical models with greater number of force 
constants [U [26J, [27] . On the base of the three-parameter simplified model for the ideal 
graphene we study further the effect of point defects on the graphene phonon spectra. 

In Section 3 we describe the contribution of isotopic defects to the graphene phonon 
DOS, analyzing separately the cases of single defect, dimer and pair of distant defects. 
In doing so we specify for graphene approaches based on the method of classic Green 
functions which were developed more than fifty years ago in [2E1 [23 EOl Ell- 
in Section 4 the same problem is considered for other point defects: substitutional 
atoms and vacancies. In these more complicated cases along with masses of defects the 
force constants between defect sites and their nearest neighbors should be changed. 

Finally in section 5 we discuss characteristic traces of considered point defects in 
optical spectra and the heat capacity of graphene. 

As illustrations we demonstrate linear in defect concentration contributions to DOS 
as more important for comparison with real experiments and estimations of defects 
manifestation in the phonon spectra. 

2. Phonons in ideal graphene 

The ideal graphene is a 2D crystal with two carbon atoms per elementary cell (further 
we will call them M" and "5") (figure [TJ. 

The equilibrium position of atoms on sublattices of the carbon plane can described 
by vectors 

R n,A = n iA a i + ™2,Aa 2 , R njB = ni,Bax + n 2 , B &2 + -^= (a x + a 2 ) 

with integer n^A, U2,a] n\,B, «2,b- The instantaneous ion configuration of the graphene 
plane is characterized by the actual positions of atoms 

R n(T = R n(T + u nCT (t), <r=(A,B), 

with time-dependent displacements 

u„ lC r(t) = M n <x,i(t)ai + u nCTi2 (t)a 2 + u na ^(t)a 3 
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around R^, where a 3 is the orthogonal to the lattice plane unit vector. Since the angle 
between translation vectors &i and a 2 is n/3, then the scalar product of displacement 
vectors u ncr and u n / CT / has the form 

(u ncr , U n / CT /) = 

Mtiff,l M nV,l + \ ( w n(T,lMnV',2 + U na t2 U n ' V ,l) + M nff ,2M n V',2 + Wncr,3 M n'o-',3- 

Therefore the kinetic energy K of lattice atoms can be written as follows 

K= ^EE^0Un, CT = 
n a 

\ E E m + «L,2 + «n ff ,l«n,,2 + «L, 3 ) , 

n cr 

where mo is the mass of the carbon atom (isotope l2 C). 

The total potential energy W of ideal graphene is modeled as a sum of three 
components determined by three different force constants: 

W = W C + W l . p . nc + W . p . nc , (1) 

where W c is the part of total potential energy depending on the relative displacements of 
neighboring atoms along the lines connecting their equilibrium positions i.e. the part of 
the potential energy due to the central forces, Wi. p . nc is the component depending only 
on magnitudes of in-plane relative displacements of interacting nearest neighbors, W , p , nc 
is the non-central component conditioned by 7r-electrons interaction, which depends on 
the out-of-plane relative displacements of neighboring atoms. For our parametrization 
the nearest neighbors to a site nA (nB) of the graphene lattice are sites (n + S)B 
((n — S)A), where S either 0- vector, or — ai, or — a 2 . Using this notation and setting 




Figure 1. The graphitic plane structure. A (full circles) and B (hollow circles) 
represent two sublattices, ai and &2 are two primitive translation vectors, |a.i.2 1 = a = 
V%b, where b — 0.142nm is the interatomic distance or the length of carbon-carbon 
cr-bond. 
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we get the following expressions for the components of potential energy in ([I]) 

W c = 4£J</l E E (Au nS ,BA, H-n+S,B - R n,A) 2 + 
n S L 

Wi /p . nc = \ J 2 E E [[Au nSiBj 4 x a 3 ] 2 + [Au nS) AB x a 3 ] 2 ] , 

T o.p.nc 



n S 



(2) 



Z J 3 E E [(Au nS ,sA ■ a 3 ) + (Au nS ,AB ■ a 3 ) ] 

n S 



with indeterminate force constants J\, J 2 , J3. 

As usual in lattice dynamics, we use further the Bloch theorem, according to which 
the atom or ion displacements on sublattice sites n for the normal modes differ only by 
phase factors exp [ik • R„], where k = k 2 ) run the Brillouin zone of the reciprocal 
sublattice. In this way we obtain that the squares of the frequencies cjj(k) for the 
different branches j of the ideal graphene phonon spectra coincide with eigenvalues of 
the dynamical matrix 



D(k) = M- 1 / 2 (JxG c + J 2 G Lp . nc + J 3 G . 
where M is the mass matrix: 
M = 



p.nc. 



rvT 1 / 2 , 



M 
M 



M = m (M i . p . + M o . p .), 



M 



1. p. 



(l \ o\ 

\ 1 







M 



o.p. 



( \ 


v° l ) 




G d c 



-M- 
2 ivj -i.p. j 



(-id 

i.p.nc 
(~i d 

^* o.p.nc 



and G c ,Gj. p . nc and G . p nc are Hermitian 2 x -block matrices of the form 



G 



with 3x3 diagonal and anti-diagonal blocks G d and G a , respectively: 

/ 1 - e ifcia 

1 - e ik2a 

3Mi. p ., Gl pMC = - (1 + e ifcia + e ifc2a ) M,. p , 
3M . P ., Gl pMC = - (1 + e ifcia + e itea ) M . p .. 

Hence for the ideal graphene plane there are six branches Wj(k) (j = 
LA,TA, ZA, LO,TO, ZO) of phonon spectra: for two of them (uza(^) and uzo(X)) 
the atom displacements are perpendicular to the lattice plane while for the other four 
branches the atoms do not come out of the plane. 



V 



o\ 
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Table 1. Exact expressions for phonon frequencies (in cm 1 ) in the r,M and K points 
of the Brillouin zone. 



u 3 



w 5 



w 6 



r 

M 
K 



' 2J 2 
m 



3(./i+2.7 2 ) 



2(./ 1 + J 2 ) 



J1+4./2 



3(Ji+2J 2 ) 
2m 



3.7i+4.7 2 



3(Ji + J 2 ) 



Setting 

F (k) = 2 [cos {kid — k 2 a) + cos k 2 a + cos fcia] , 

Fx(k) = 12(J2 + 2Ji J 2 + 2J|) + F (k)(^ + 8JiJ 2 + 8J|), 

= J x 2 + 16 JtJ 2 + 16 J|, X 2 = J\ - 8Jl J 2 - 8J|, (3) 

F 2 (k) = {18 + 2Xi + 4X 2 cos (k±a — k 2 a) (cos k±a + cos A; 2 a) + 

+4 cos kid cos fc 2 a [X 2 + X\ cos (&ia — k 2 a)] — 6 J 2 F (k)} 2 , 

we get the following expressions for the eigenfrequencies Wj(k): 

1/2 



(k) 



3± V / 3 + F (k) 



WLA,TA,LO,To(k) 



3(Ji+2J 2 ) 
2m 



±^ v /F 1 (k)± v / 2JiF 2 (k) 



- q l/2 



(4) 



The last expressions take a very simple form at the high symmetry points T (ka = (0, 0)), 
K (ka = (|7r, |7r)) and M (ka = (tt, tt)) of the graphene first Brillouin zone (Table [l]). 

In the vicinity of T-point there are three acoustic branches of the graphene phonon 
spectra, the eigenfrequencies of which in the limit k — > do not depend on the 
direction of propagation. They can be distinguished by polarization of oscillations as 
the longitudinal branch u>la = cla^ + 0{k 2 ) with longitudinal sound velocity 



(Ji + J 2 )(Ji+4J 2 ) 

cla = a* ; (5) 

y 6(Ji + 2J 2 )m 

the transverse in-plane branch uta = CTAk + 0(k 2 ) with in-plane transverse sound 
velocity 



cta 



^J X J 2 +4:Jl 



6(Ji + 2J 2 )m ' 

the transverse out-of-plane branch ujza 
sound velocity 



cza 



3mn 



(6) 



czAk + 0(k 2 ) with transverse out-of-plane 



(7) 



By definition the density of states p(co) of any oscillatory system is given by 
expressions 

p(u) = 2uTZ(u 2 ), tu>0, 



W(A) = E*(A-wJ) = J%E 



e4,o 



(A-^) 2 + £ 2' 



(8) 
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Table 2. Theoretical phonon frequencies for the T, M and K points of the Brillouin 
zone. 

Phonon frequencies (cm -1 ) 

UJl UJ 2 W 3 CJ 4 W 5 UJq 

r 840 1620 
M 485 686 828 1031 1249 1392 
K 594 1014 1146 1263 



where uj u are eigenfrequencies of the system counted with respect to their multiplicities. 
In our case the normalized to 6 ( = the number of degrees of freedom per elementary 
cell) DOS is given by the formula 

■k J a it I a g 

= - iim i^-) 2 r r J2- £ —— 2 — dk 1 dk 2 = 

7r/a n/a 

Mimlm(^) 2 / / Tr[D(k)-(A + k)J] _1 (iMA;2, ( 9 ) 

e ^ —ir/a ~n/a 

oo 

/ H(\)d\ = 6 



where X is the 6x6 unity matrix. 

In order to obtain the phonon dispersion curves and DOS for ideal graphene on 
the base of expressions Q and (12) concrete numerical values of the force constants 
Ji,J%, J3 are necessary. We derived these constants by the least-squares method using 
experimental values of graphene phonon frequencies at the T, K and M points measured 
using the method of inelastic x-ray scattering. From now on we assume that 

Ji = 3.79 x 1(T 21 J/m 2 , 

J 2 = 6.89 x 10- 21 J/m 2 , (10) 
J 3 = 2.36 x 1Q- 21 J/m 2 . 



With the constants (10) we get the values of phonon frequencies at the T, K and 



M points of Brillouin zone (Table |2j) and the sound velocities (which are in a good 
agreement with [U El [6] ) 

cla = 18.4 km/sec, 
cta = 16.5 km/ sec, 
cza = 9.2 km/ sec. 

For demonstration of how the fitted three-parameter model works we give below the 
phonon dispersion curves (figure [2]) and density of states (figure [3]) for ideal graphene. 
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3. Isotopic defects 

The kinetic energy K and potential energy W of the non-ideal graphene can be written 
as the following quadratic forms 

1 a 
1 cr 

- B 3 

W ({11,(1)}) = S Civw^'i'O'KiG), 

l',I cr',cr=Aj',j = l 

where la enumerate the graphene lattice sites. From now on we will denote by 9JT and 
£ the matrices of the quadratic forms 2K and 2W, respectively. We will call 9JT the 
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Figure 2. Theoretical phonon frequency dispersion curves of the graphene along 
r — M — K — r directions and experimental data for graphite (solid circles) [4]. LA 
(LO), TA (TO) and ZA (ZO) are longitudinal, transversal and out-of-plane acoustical 
(optical) branches respectively. 




Figure 3. The phonon density of states of the infinite graphene plane without defects. 
The phonon frequencies at the high-symmetry points of the Brillouin zone are labeled 
by r„M, and K,. 
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mass matrix and £ the force constant matrix. For the case of the ideal graphene SOT is 
the block diagonal matrix S0T with equal diagonal blocks 

Ml M = M (11) 

and €. is the matrix of the doubled sum of the quadratic forms ([2]). For graphene with 
substitutional defects the mass matrix is block diagonal and all its diagonal blocks have 
form (11) but with masses of substitutional ions instead of mo for some la. 

Let SOT be the mass matrix of a non-ideal graphene plane with some host 12 C atoms 
replaced by other stable carbon isotopes. Note that the block-diagonal matrices SOTo and 
SOT commute. Let po(u) and Pvk{lo) denote DOS of carbon planes with mass matrices 
SOTo and SOT, respectively. As to 2 in ^ are eigenvalues of the matrix 

D =50H£9jri 

or the linear pencil € — zSOT, i.e. u 2 are those values of z, for which £. — z$R is non- 
invertible, then 

n(u 2 ) = 1 lim ImTr [D - (to 2 + ie) 3}~ l = 

W ei ° (12) 
Mim ImTr[£- (w 2 + i£)50T]- 1 50T, V ; 

where 3 is the unity matrix. 

Note that for any positive definite matrices € and SOT of any finite order and non-real 
z we can write 

Tr ([£- zSOT^SOT) = ~ Tr In (C - zSOT) = ~ lndet(£ - zWl). (13) 
Therefore from ([8| and (p~3l we get 



or 



PanM = p (u) - — lim Im{-^lndet [(€ - zM)(€ - ^COTo)" 1 ] 1 (14) 



Pot(w) - po(w) = 
-f limlm£ lndet [3 — z (SOT — 9JT ) (£ - ^SJTo)" 1 ] \ zMl 

limlm£ lndet [3 - z (SOTSOTq 1 - 3) (2) - z^ 1 ] | Z=U)2 _ 



(15) 



From now on we will assume, for simplicity, that only one species of isotopes with 
mass m can replace the host atoms of the ideal graphene lattice and assign to each site 
la of non-ideal graphene the occupation number 

{0 if the host carbon atom seats at ler; 
1 if the carbon isotope is there. 

We can express diagonal blocks of SOT in terms of n\ a as follows 



M la . 



la 



m-m 

1 H nia- 

m 



M . (17) 



By (17) (SOTSOTq 1 — 3) is a block diagonal matrix with 3x3 diagonal blocks of the form 

m — mo 

Ma,la = P ■ ni a l, pL= , (18) 

m 
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where I is the 3x3 unity matrix. 

In view of the obvious property of occupation numbers n\ a = n\ a we can represent 
Pb(w) formally as follows: 



PoM + yi Y,-i( UJ '^ a ) n ^ + | E H 2 (u;;la,lV)n lo .7ii/ CT /+ 

1(7 (19) 
± E S 3 (o;;la,lV,l / V)ni (r ni' CT /ni// CT // + ... 



In fact, the decomposition (19) is an identity, which holds for any number of isotopes (or 



other point defects) and their distribution over the graphene lattice sites. Particularly, 
if there is only one isotopic defect located at the lattice site locr , that is if n\„ = 5kr,i 0(T0 , 



then according to ( 19 ) 



Si (a;; l a ) = l a ) - Po{u), (20) 

where p(cj; loa ) is the DOS of graphene lattice with a single defect at the site l a . In 
much the same way we find that 

2 2 (w; la, IV) = p(u; la, IV) - Ex(u; la) - Ei(u>; lV) - p (w), (21) 

where p(u; la, IV) is the DOS of graphene lattice with only two isotopic defects at the 
sites la and IV and so on. 

Let T(z) denote the 3x3 diagonal block of (D — zJ)' 1 with some index la. By 
the translational and point symmetry of graphene plain T(z) does not depend on la. 



It follows from (15) and (17) that the coefficients Hx(a;;loao) in (19) actually do not 



depend on 1 and a, 



(w;Vo)) 



lim Im— In det [I — uzT(z)] 

71 elO dz 



Note that for 3 x 3 blocks d W /„/ of 2) = [m o 2 £9Jt 2 J we have 

dl CT; lV' — di-i^oo-'; dn/o-.oo- = di_iv ;0 T, o ^ r. 

Making use of the explicit expression for the dynamical matrix D(k) = {T> (TT {k)) 2 (7T=l of 
the ideal graphene plane 



(k) = -Vd Mr e- M , 
m 



where k is the wave vector from the Brillouin zone, we obtain that 

k 

where N is the number of unit cells in the periodicity area. Hence 



- / \ 2u d 
'Moj) = lim lm — In det 

7T £+0 dz 



l -/*4£(p>( k )-* i r 1 ) t 



z=ui^+ie 



Proceeding in the same fashion and setting 

r -(^ 1 ) = ^E eik ' 1 ([ I? ( k )- 2j r 1 ) c 



(22) 



(23) 



(24) 



(25) 
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we find that 
~ 2 (cu; kr, IV) = S 2 (w; 1 - IV, Or) = 

- 2Sl H - ^limlm-f Indet f n "^^J 7 1 ■ (26) 

V ; n eio dz \ -zfiT ra {z;V -\) I-znF(z) J 




Figure 4. The influence of 2% aluminium defect concentration (isotopic case) on the 
ideal graphene phonon density of states. The dashed, solid and dotted lines represent 
DOSs of the ideal graphene, defect graphene (2% of Al atoms) and defect contribution 
respectively. 




Figure 5. The influence of 20% nitrogen doping (in the case of dimers) on the ideal 
graphene phonon density of states. The dashed, solid and dotted lines represent DOSs 
of the ideal graphene, defect graphene (20% of N atoms) and defect contribution 
respectively. Two gaps in the density of states near 5 and lom 6 are labeled by 
the arrows. 



Phonons in graphene with point defects 12 
We pay now attention to the fact that for the graphene plane of finite size 



p m (oj), Pq{uj) in (19) have order iV while the coefficients Hi (a;), H 2 (a;; la, Or), ... are 
finite quantities. Let us consider the equilibrium distribution of defects (isotopes) with 
concentrations c and the pair distribution function of defects 

g ar (l) = lim N(n la n 0T ), 

V— s>oo 

where angle brackets denote thermal (or other) average. Then for the phonon spectral 
densities per unit cell of the graphene plane 

A*(w) = lim -{p A (u)),St=fm,0, 

TV— s>oo iV 

by @, ([24]), (|26|) we have 

A m (u) = A (w) + c2i(w) + ^c 2 ^ ~ 2 (cu; la, Or)g ar (\) + ... (27) 

1(7T 

If for some reason there is a visible concentration p of isotopic dimers on the graphene 
plane, that is cells occupied by pair of isotopes, then their contribution to A«jk(u;) is 
(figure [5} 

A«ut(w) = ... + p • H dim (cu)..., 

I-^r(^) -z/iL^^O) \ (28) 



-dim { U ) 



^limlm-f lndet . 
77 eio dz \ -zpT Ta (z;0) I-zpT(z) 



4. Substitutional atoms 

Let us consider now a non-ideal graphene plane with low concentration cx of impurity 
atoms X replacing carbon atoms at some lattice sites. Applying the same arguments 
as above we can again assert that in this case the phonon spectral density per unit cell 
Ax(u) can be written as follows: 

A X H = AqM + £*6i(o;) + ^4 Yl Q ^ la ' 0r )3-(l) + o(c 2 x ), (29) 

\<JT 

Here 

QiM = p x {w) - po{u), (30) 

px(w) being the DOS of graphene lattice with a single defect at some site lo^o; 

Q 2 (cu; la, Oa') = Px (oj; la, Oa') - 2Q 1 (u) + p Q (u), (31) 

Px{u; la, Oa') being the DOS of graphene lattice with only two impurity atoms X at 
some sites Ya, I" a' such that 1 = 1' — 1". Denote by VJtx the mass matrix and by £x the 
force constants matrix of the graphene lattice with a small number of impurity atoms 
X replacing carbon at some sites and put 
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(32) 



As a starting point for calculation of Qi(u), 02(w; la, Oa'), ... we make use of the 
relation 

Px(uj) - PoM = 
- 2 f Hmlm£ lndet [3 + {5£ x - z6Wl x ) (£ - zWlo)' 1 ] \ z=u)2+iE . 

Note that as well as for the graphene lattice with only isotopic defects the matrix 

T x (z) : = 9JH (5£ x - z5Wl x ) 3jH 

can be represented in the block form, elements of which are 3x3 matrices being 
enumerated by pairs of the multi-indices la, IV that enumerate the graphene lattice 
sites. Actually, only those blocks of T x (z) are non-zero, for which each of the indices 
la, IV belongs to the subset of indices D of those enumerating either sites occupied by X- 
atoms or their nearest neighbors. Therefore if there are N x impurities in the periodicity 
area, then the rank of T x (z) does not exceed 12N X . For these reasons the determinant 
in (32) can be replaced by the minor determinant obtained from det(3 + ...) by deleting 
all rows and columns with "numbers" other then those from D. By this argument 
setting 

Q(z) = [T aT (z; I - r)L,i/ re 3 , T X (z) = [T X {z)la,VT\ a>VT ^ 

we can write 

p x (u) - p (u) = 
-f hmlm £ lndet [t x {z) \T x {z)^ + Q(z))- 1 ' 

Let us assume for certainty that the single defect is located at the site OA and 
therefore its "normal" nearest neighbors are located at the sites li-B, hB, 1%B, where 
li = 0, 1 2 = —a!, 1 3 = — a 2 , respectively. In this case the subset introduced above 

D = {0A, hB, l 2 B, l 3 B, hB}. 

We may count further for brevity instead of OA and s instead of l s B, s = 1,2, 3. 

With this enumeration the 12 x 12 matrices T x {z) and T x {z)~ l can be now 
represented as 4 x 4 block matrices 



(33) 



z=u) 2 +ie 



( 



Tx(z) 



Tx{z) 



-1 



-z5G + $G S 

s=l 

-5Gi 
-5G 2 
-5G 3 

/ I I I I \ 

! I I I I 
^ I I I I 
I I I I 



-5Gi 6G2 



6G 1 





V 



+ 





6G 2 


/ 
5G^ 

y 



5G 3 * 




5G 3 ) 
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where 



5M = pi 



(35) 
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Figure 6. The influence of 2% aluminium defect concentration (isobaric case, 
Jj = O.SJ^o) on the ideal graphene phonon density of states. The dashed, solid and 
dotted lines represent DOSs of the ideal graphene, defect graphene (2% of Al atoms) 
and defect contribution respectively. 



p(w), cm 




Figure 7. The influence of 2.5% of vacancies on the ideal graphene phonon density 
of states. The dashed, solid and dotted lines represent DOSs of the ideal graphene, 
defect graphene (2.5% of vacancies) and defect contribution respectively 
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By (34) and (35) we have 



detTx(z) = - (detppG s 



s=l 



Substitution of the above given expressions (|34j) - (36) and 

/ T AA (z;0) T AB (z;0) T^^a^ 
T BA (z;0) T BB (z;0) T BB (z;&i) 
r BA (z;-a.i) r BB (^;-ai) T BB (z;0) 
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T BB (z;0) ) 



into (33) gives the sought expression for the "coefficient" @i(ui) in (29). 



In the partial case of "isobaric" defect (/i = 0) we have 
e 1 (u) = -^limlm£lndet @ + T Q(r" 1 
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(38) 



(39) 






5G 3 ) 

Not much remains to add to prove that for the vacancies at the graphene lattice 
sites the corresponding coefficient Q\[bS) is given by the expression 
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(40) 



where bG\,bG\,bG\ are the particular cases of expressions (36) with AJi,AJ 2 ,AJ 3 
replaced by —J\,—Ji,—J%, respectively 



5. Discussion 



The proposed model for the description of phonon spectra in graphene, where the only 
nearest neighbors interaction is accounted for, gives rather simple explicit expressions 
for the phonon dispersion curves and reflects quite satisfactory their main features for 
graphene and graphite in the whole frequency range. Certainly, containing only three 
force constant such a model can't reproduce all the specific features of graphite phonon 
spectra. For example, it cannot in principle reproduce the low-frequency bending mode 
with bj ~ k 2 and the overbending of in-plane optical modes near the T point so the more 
that the latter is assumed to be caused by the electron-phonon interaction [201 E3J El] ■ 
Nevertheless, the three force constants of the model being chosen as fitting 
parameters to obtain a good coincidence with known values (inelastic x-ray scattering 
on graphite |4j) of optical modes at the symmetrical points of the first Brillouin zone 
r 5 6 , K 6 , K4 5, K 3 , M 6 , M 5 (for in-plane modes) and r 4 , Mi (for out-of-plane modes) give 
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satisfactory (and sometimes rather good) quantitative agreement with experimental data 
for graphite throughout the whole frequency range. In particular, the value of in-plane 
transverse sound velocity agrees closely with the experimental value and results of the 
more sophisticated theoretical models |H El El- We took slightly overstated in-plane 
T-point frequency wr 56 = 1620cm -1 (instead of a value in between 1580 - 1595 cm -1 
according to known experiments) for better fitting of experimental data in a larger part 
of Brillouin zone and for a formal account of the in-plane optical modes overbending 
near the T-point. Since the experimental data for K and M points of graphene phonon 
spectrum are rather poor and the selection rules are strictly obeyed due to the high 
crystalline quality of graphene we take the graphite phonon spectrum as a target of the 
fitting procedure, more so that the phonon dispersion curves for graphite and graphene 
are similar in a large part of the Brillouin zone [23] . 

Another visible discrepancy between our theory and experimental data (in addition 
to those mentioned above) is the sufficiently lower in frequency in-plane longitudinal 
acoustic (LA) branch (figure [2]) in the T — M direction and as a result its incorrect 
trend between M and K points. The analytically obtained frequencies uk 6 and um 4 are 
connected in the model by the relation Uk 6 = y§ % which is different from the x-ray 
experimental result uk 6 ~ [!]• However, this relation is in a better agreement with 
Raman experiments which give ujk 6 ~ V r b6wM 4 [32] - Also for the o;r5 i6 and ujk 4S we get 
cor 56 = v2Wit4 S ) while x-ray and Raman experiments give 1.32 and 1.25, respectively, 
(instead of 1.41 (or 1.38 with account of our Co>r 56 overestimation)). The next interesting 
observation is that the com & , ^Ki B , ^k 6 > ^rs.e ratios to % 4 which are drastically different 
from the x-ray experiment results are in a good agreement with Raman experiments 
despite of essential difference in absolute values of some of those frequencies [32J. It 
is worth mentioning that the similar relations and ratios were obtained within a more 
detailed model with up-to-third nearest neighbors interactions accounted for [26] . 

Note that since the most substantial spectrum changes due to defects were expected 
near the Van Hove singularities, all of which are reproduced at the proper places by our 
simple model, it is evident that further improvements and refinements were not dictated 
by the objectives and needs of this work. 

The study of the influence of point defects on the graphene phonon spectra on the 
basis of our simplified model for the ideal graphene lattice oscillations showed, as might 
be expected, that the isotopic defects slightly downshift Van Hove singularities in the 
phonon DOS for heavier than carbon atom defects and upshift them for light defects. 
For the defects with mass of 11 u.amu (like n B) the additional Van Hove peak appears 
on the upper edge of the ideal graphene phonon spectra. For lighter defects this peak 
splits out and corresponds to the localized oscillation mode with frequency above Co>r 56 - 
Such peaks may be related to the so called D' (u = 1620 cm -1 ) and 2D' (u = 3250 cm' 1 ) 
bands of the graphite Raman spectra [331 El] and were directly observed for 2.66% boron 
doped graphene [35] . 

Substitutional defects, for which the force constants for defect-carbon bonds are 
weaker than that for carbon-carbon bonds, may drastically change the phonon DOS in 
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the low-frequency region (from 60 to 200 cm -1 ) depending on the defect atom mass and 
weakened force constants. Note that a trace of an oscillation mode with a frequency of 
about 100— 120 cm -1 near the T-point was observed in graphite pfl |25j [36]. It turns out 
in some cases that a defect concentration even of several per cent may result in 100% 
low-frequency DOS increase in intervals wider than 100 cm -1 . Obviously, the changes 
of the phonon low-frequency DOS stipulated by defects should be manifested in the low 
temperature specific heat of non-ideal graphene. The isotopic defects with greater then 
12 C mass lead to increase of the specific heat in proportion to the defect concentration 
(and temperature). Light defects slightly decrease the DOS and low temperature specific 
heat, but the difference from ideal graphene is sufficiently less than in the case of heavy 
defects. For substitutional defects the dierence in density of states in the low-frequency 
region between graphene with defects and ideal graphene is much higher then in the 
case of isotopic defects and leads to significant changes in the specific heat (of order of 
10%). The specific heat of 3% Al doped graphene, AJj/Jj = —0.25, differs by more 
than 10% of that of ideal graphene. 

Apart from slight shifts of Van Hove points, significant changes in the DOSs are 
observed for any types of point defects near Uk 12 and % 46 points where oscillation 
modes are double-degenerated. For heavier isotopic defects the DOSs have more sharp 
peaks at the mentioned points. For substitutional atoms the changes near the u)k 12 
and ujk 4<5 are similar to those for purely isotopic defects (figure [6]). For example, the 
deviation from the ideal DOS can reach as high as 30% for 2% Al doped graphene 
(figure [4]). In comparison with isotopic defects the same concentrations of vacant lattice 
sites lead to a greater DOS rise for frequencies from uk^ 5 to com s and to sufficient DOS 
decrease in the intervals (% 2 ,ur 4 ) and (ujm 6 , w r 5 , 6 )- That lowering may be of the order 
of several per cent for the wide frequency range (figure [7]). Both types of defect also 
lead to the high additional peaks in the (see also [37]) and ujk 4:5 (figures [I| [6j). 

The singularities near the 0Jx ia and ujk 4 5 may indicate the splitting of corresponding 
modes in the vicinity of the K-point and in the case of high defect concentration the 
appearance of gaps in the density of states near cok 45 and % 6 is evident (figure [5]). 

The vacant lattice sites except for the singularity at the ujk A5 point cause an 
additional resonant states peak with frequency between 0Jx ia and ojm 2 (figure [7]). 

In the ideal graphene the momentum conservation allows only one single-phonon 
Raman process which corresponds to the emission of the optical phonon with zero wave 
vector and frequency near 1580cm- 1 (u;r 56 ) [38J . For non-ideal graphene there is also the 
so-called double-resonant presumably defect-induced peak (the D band) at 1350cm -1 
corresponding to emission of an optical phonons with wave vector near K points of the 
Brillouin zone j38j [39J, HQ]. Within the commonly accepted interpretation, impurities 
only assist the photons scattering on single intervalley phonons. On the other hand, 
our analysis shows that all considered types of point defect cause additional Van Hove 
singularities of the DOS in the neighborhoods of the K and M points of optical branches, 
and that local oscillation modes may induce additional peaks of the Raman spectra. 
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